\documentclass[orivec]{llncs}

\usepackage{headers}
\usepackage{macro2}
\usepackage{xspace}
\usepackage[noend]{algpseudocode}
\usepackage{algorithm}


%opening

\title{Activity Networks with Delays}
\subtitle{application to toxicity analysis}

\author{%
 Franck Delaplace\inst{2} \and
  Cinzia Di Giusto\inst{1} \and
  Jean-Louis Giavitto\inst{3} \and
  Hanna Klaudel\inst{2} \and
  Antoine Spicher\inst{4}
}

\institute{%%
  I3S, Université de Nice Sophia-Antipolis, Nice, France \and
  IBISC, Université d'Evry - Val d'Essonne, Evry, France %
  \and %
  CNRS, IRCAM, UPMC - UMR 9912 STMS \& INRIA MuTAnt team, Paris, France 
  % 
  \and %
  LACL, Université Paris-Est, Créteil, France}


\begin{document}
%\widowpenalty=0
%\clubpenalty=0

\maketitle 
\begin{abstract}
      \andy, Activity Networks with Delays, is a discrete time framework aimed to the qualitative modelling of time-dependent activities such as biological  networks and regulatory pathways. 
     
      Activities involve entities playing the role of
      activators, inhibitors or products of biochemical network operation. 
      Activities may have given duration, \ie the time required to obtain results. 
      An entity  may represent an object (\eg an agent, a biochemical species or a family of thereof) with a local attribute, a state denoting  its level (\eg concentration, strength). 
      Entities levels may change as a result of an activity  or may decrease gradually as time passes by. 
      
      The semantics of \andy is formally given via  high-level Petri nets ensuring this way some modularity.
      As main results we show that \andy systems have finite state representations even for potentially infinite processes and that  \andy's concept of time has a direct counterpart in timed automata. 
      
      These results together with a modular and concise syntax  make \andy suitable for an easy and natural modelling of  time-dependent (biological) systems. 
      We conclude our paper with a general discussion on toxic behaviours, in particular we present a classification of toxicity properties and give some hints on how they can be verified with existing tools on \andy systems. 
\end{abstract}
 

\section{Introduction}

Activities are a central concern in many domains such as business (workflows, supply chains), transaction modelling (database transactions, web services) or biological systems (biological networks, regulatory pathways).
%, to mention some. 
These activity-driven domains, although quite different, share foundational principles and concern: the causal relationships among activities. It deals with questions such as,  how long does it take to complete an activity, when it can be performed,  
 or whether it can be delayed or even ignored. 
%Control-data, instead,  addresses how data can play a role in promoting or  preventing an activity.\looseness=-1

Our main motivation comes from biological systems where data are species and the activities  the reactions among them but we believe that our model could be adapted also to handle other  domains as the ones mentioned above. 
Our proposal stems from reaction systems %as defined in 
\cite{DBLP:journals/ijfcs/BrijderEMR11}, a formalism based on  \emph{reactions},
each defined as a triple $(R, I, P)$ with $R$  set of
\emph{reactants}, $I$  set of \emph{inhibitors} and $P$  set of
\emph{products}, and $R,I$ and $P$  taken from a common set of
\emph{species}. Reaction systems are based on three basic assumptions:
\begin{enumerate}[(i)]
 \item \label{item:react} 
 a reaction can take place only if all the reactants involved are available but none of the inhibitors are; 
 \item \label{item:quantity}  if a species is available then a sufficient amount of it is necessary to trigger a reaction; 
 \item  species are not persistent: they become unavailable if they are not sustained by a reaction.
\end{enumerate}
%
We complement this model by adding several important features.
We allow a richer description of species states using attributes having potentially multiple values  and introduce timing aspects.  
More precisely, the resulting formalism, \emph{Activity Networks with Delays} (\andy), 
targets systems 
composed of various \emph{entities} (\ie species in reaction systems) that evolve by means of time-dependent activities.
Each entity is characterised by one attribute,
called \emph{level}, representing rates, activation/inhibition status or expression degrees
(\eg low, medium, high). 
\emph{Activities} are rules describing the evolution of systems, involving (as for reaction systems) three sets of entities:
\emph{activators}, which must be present, 
\emph{inhibitors}, which must be absent,  
and \emph{results} whose expression levels have to be modified (increased or decreased).
The introduction of time concerns both activities and entities. Activities have a duration and entities are subject
to a \emph{decay} process or aging that models the action of a
non-specified environment: \ie levels decrease with time
progression.
Another difference is in the semantics of activities: while in reaction systems a maximal concurrency model is considered, here we  
adopt two types of activities:  \emph{mandatory} and \emph{potential} ones. The former set of activities, once enabled, must take place altogether in the same time unit (in maximal concurrency) while the latter, non deterministically one at the time, may be performed or not.  
Following \cite{Wang04formalverification}, as in our modelling all entities have discrete states and share the same global clock, 
it is reasonable to work with discrete time constraints.\looseness=-1 


 Our main objective is to provide effective tools for the
 description and  the understanding of the mechanisms and of the
 structural properties underpinning biological interactions networks. In
 biology, ordinary differential equations (ODE) remain the
 predominant modelling methodology\cite{alon2006introduction}. Such models present however some
 drawbacks: they need a precise quantification of parameters, 
random interactions are difficult to model outside an averaged
 approach, %
 system openness and non-linearities make hopeless the existence of
 analytic solutions, %
 and numerical methods -- often the only sensible option -- are mainly
 descriptive and cannot form the basis for an explanatory theory. 
 On the other hand, many
 languages and tools  have been developed  to build and explore discrete qualitative
 models in various application domains. Such models are rough
 abstractions of real world processes but they make possible to unravel
 the entangled causal relationships between system's entities.
 %%
 More specifically in systems biology, formalisms like Petri
 nets~\cite{DBLP:journals/nc/BaldanCMS10}, Boolean networks~\cite{thomas1973boolean},
 process algebras~\cite{journals/tcsb/Cardelli05} 
or rewriting
 systems~\cite{giavitto04a} (to cite a few) have been used with
 undeniable success to understand the causal links between structure
 and behaviour in molecular networks.
  Nonetheless, the management of time is somewhat problematic in all of
 the previous approaches. 
 In general, timing aspects are either disregarded or handled at a primitive level,
 which leads to expressiveness problems for the modeller. 
 For example, modelling a synchronous evolution directly in Petri nets
 requires a coding that is not necessarily natural for a biologist.
 % In general, time aspects are either disregarded or they are handled with stochastic or continuous methodologies. 
 This is why here we propose  a qualitative discrete formalism for biological systems that %explicitly takes into 
 provides, in particular,  a direct account of timing aspects natural in biology such as activity duration and decay. 

We complement our designing process with a general discussion on toxicity properties and how they can be tested against \andy systems. In this respect, we aim at providing a methodology rather than a precise technique of testing.

\paragraph{\bf Organisation of the paper.}  This paper is the extended and revised version of BIOPPN workshop paper \cite{CDGHKFD14}. The main difference with respect to \cite{CDGHKFD14} is the introduction of time aspects in \andy systems and a more mature taxonomy of toxicity properties.

%The paper is organized as follows: 
%%
Section \ref{sec:reaction} introduces the principles behind \andy
networks
and gives a formal definition of \andy in terms of high-level Petri nets. It also states the main results of the paper  addressing the finiteness of state space and the expressiveness of the formalism in terms of timed automata.
%%
The encoding of an \andy network into timed automata is exponential 
showing the expressive relevance and the concision of the original
\andy's notation.
%%
In Section \ref{sec:toxic} we discuss an axiomatisation of toxicity properties and finally Sections \ref{sec:related} and \ref{sec:concl} discuss related works and conclude.
%with some perspectives opened by the first simulation and model-checking
%tools developed for \andy.

An Appendix with additional material: the details of the encoding into timed automata and a supplementary  example on the assimilation of aspartame in human body is added for the reviewer convenience. 




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{Activity Networks with Delays}\label{sec:reaction}

In this section we give the syntax and   
semantics of \emph{Activity Networks with Delays} (\andy) which is composed of a set of entities
\entities driven by a set of activities. We assume to work with a discrete notion of time represented by a tick event.


\paragraph{\bf Entities.} 
%%
Each entity, $\entity \in \entities$, is associated to a finite number of levels $\setlev_{\entity}$ (from 0 to $\setlev_{\entity}-1$), that 
in general represent ordered expression degrees (\eg low, medium,
high) and refer to a change in the entity capability of action. A decay duration is associated to each level of an entity $e$ 
through the function $\life_\entity:[0 \mydots \setlev_\entity-1] \to \nat^+ \cup \{\omega \}$ 
allowing different duration's for different levels and where $\omega$ means unbounded duration. %Decay is a  natural number or $\omega$ (unbounded),
If $\life_\entity(i) = \omega$, then level $i$ is \emph{permanent} and it can be modified
only by an activity. Otherwise, if $\life_\entity(i) \neq \omega$, the presence of the entity
at level $i$ is \emph{transient}: once entity $\entity$ at level $i$ has passed $\life_\entity(i)$ time units its level will pass to $i-1$, \ie it decays or ages gradually
as time passes by.
We fix $\life_{\entity}(0) = \omega$ as no negative level is allowed. 
 In principle, this substitutes a set of unspecified activities accounting for the action of an underspecified   
environment that consumes entities.  
A state $\state$ of an \andy network assigns to each entity $\entity \in \entities$ a level $\eta \in [0..\setlev_{\entity}-1]$. The initial state $\state_0$ sets each entity $\entity$ to a given level $\level{\entity}$. 
 

\paragraph{\bf Activities.}
The evolution of entities is driven by \emph{activities} 
 $$ \rho ::= \activ{A_{\rho}}{I_{\rho}}{\dur{\rho}}{R_{\rho}}$$  
where $\dur{\rho} \in \nat$ is the activity duration, $A_{\rho}$ (activators) and $I_{\rho}$
(inhibitors) are sets of pairs $(\entity, \level{\entity})$, and $R_{\rho}$
(results) is a non empty set of pairs $(\entity, v)$, with $\entity\in
\entities$, $\level{\entity} \in [0 \mydots \setlev_\entity-1]$ and $v \in
\mathbb{Z}$. 
Entities can appear at most once in each set $A_{\rho}$, $I_{\rho}$ and $R_{\rho}$, cf., Remark \ref{rem:def} below. 
We write $e \in A_{\rho}$ to denote  $(e, \cdot) \in A_{\rho}$, similarly for $I_{\rho}$ and $R_{\rho}$. We omit index $\rho$ if it is clear from the context.

An activity $\rho$  can take place only if it is \emph{enabled}, this depends on
 the activator and inhibitor
levels appearing in  $A_{\rho}$ and $I_{\rho}$. Nonetheless, observing only the current level of each entity is not enough to
decide if an activity is enabled: such approach will miss the handling
of decays and duration. 
%%
To take into account a notion of elapsed time, an activity $\rho$ is
enabled (\ie may happen) if and only if:
\begin{itemize}
 \item for each entity $\entity_a$ of the activators set $A_{\rho}$ (\ie $(\entity_a, \level{a}) \in A_{\rho}$), $\entity_a$ is available at least at level $\level{a}$ for
  the whole duration $\dur{\rho}$;
 \item for each entity $\entity_i$ of the inhibitors set $I_{\rho}$ (\ie $(\entity_i, \level{i}) \in I_{\rho}$), $\entity_i$ is available at a level strictly inferior to $\level{i}$ for
  the whole duration $\dur{\rho}$;
\end{itemize}
If an entity is an activator (resp. an inhibitor) at level
      $\ell$ for some activity $\rho$, it is also an activator (resp. an
      inhibitor) for $\rho$ at all levels $l \ge \ell$ (resp. $l \le \ell$).
When an activity takes place, entities' levels appearing in
the results are updated with the increment or decrement (up to the range of entity levels) $v$ specified in the rule. Duration $\dur{\rho}$ characterises the  number of ticks required for yielding change of
levels of results, $\dur{\rho} =0 $ denotes an instantaneous activity.
%%


\begin{remark} \label{rem:def}
Notice that:
\begin{itemize}
 \item  An entity may be both in $A$ and $R$, \eg as in an
      auto-catalytic production where reaction products themselves are activators for the reaction. Similarly, an
      entity may be both in $I$ and $R$: thus representing self-repressing activities.
      
 \item Activators and inhibitors are not ``consumed''.

      
\item An entity can appear in the same activity simultaneously as an
      activator and as an inhibitor but we require them to
      occur with different levels. For instance, the rule
      $$
      \activ{\{({\entity},\level{})\} \cup A}{\{({\entity}, \level{}')\} \cup I}{\dur{}}{R}
      \qquad  \text{where~}   \level{} < \level{}'
      $$
      requires that the activity takes place only if the level
      $\lev_{\entity}$ of ${\entity}$ belongs to the interval $\level{}
      \leq \lev_{\entity} < \level{}'$. In particular, if ${\entity}$
      has to be present in an activity exactly at level $\level{}$,
      ${\entity}$ should appear as an activator at level $\level{}$ and
      as inhibitor at level $\level{}' = \level{}+1$.

      
\item The set of activators and inhibitors $A \cup I $ is allowed to be
      empty. In this case, the activity is constantly enabled. This
      accounts for modelling an environment that  continuously
      sustains the production of an entity. In this case a duration
      $\dur{} \not= 0$ may be used to account for a start-up time.

\item Once triggered an activity $\rho$ must wait at least $\dur{\rho}$ before being enabled again\footnote{This is a non restrictive hypothesis and by slightly modifying Definition \ref{def:rs} one can obtain different semantics that could depend on the specific modelled scenario: \eg once enabled an activity $\rho$ can be triggered an unbounded number of times or  once per time unit.}.
      
\end{itemize}
 
\end{remark}

It is important to note that when an activity is enabled, it is not
necessarily performed. The idea is to make a distinction between
two types of activities: 
\begin{enumerate}
\item \emph{\potential} activities, denoted $\alpha$ in a set $\Reac$:
      if enabled, they \emph{may} occur in an interleaved way thus the result
      depends on the order they are performed. It may happen that no \potential activities occur even if they are enabled;
\item \emph{\obligatory} activities, denoted $\beta$ in a set $\Syn$: if
      enabled  they \emph{must} occur within the same time unit. All enabled
      \obligatory activities occur simultaneously.
\end{enumerate}

As all \obligatory activities $\beta_i$ are triggered
simultaneously, all updates related to an entity $\entity$ are
collected and summed up before being applied to the current level of
$\entity$. 

\begin{definition}[\andy network]
An \emph{\andy network} is a triple $(\entities, \Syn, \Reac)$ where
$\entities$ is the set of entities, $\Syn$ is the set of \obligatory
activities and $\Reac$ is the set of \potential ones.  
\end{definition}

An {\andy network} evolves by triggering enabled activities in two separate phases implementing one evolution step:
\begin{enumerate}[Ph. 1]

\item \label{item:reaction} %%
      In this phase (between two ticks) \textbf{\potential activities} that are
      enabled may be triggered. Their action occurs non-deterministically in an interleaved
      way. Any potential activity is triggered at most once in this phase. 

\item \label{item:clock} %%
      This phase is the tick transition: the effect is that all  concerned entities must decay and  all enabled \textbf{\obligatory
        activities} are performed simultaneously.
      %% 
      Notice that as activators and inhibitors are not ``consumed", there
      are no conflicts and so there is a unique way (per time unit) to
      execute all enabled \obligatory activities.
\end{enumerate}



\begin{example}[Repressilator]\label{ex:rep}
The repressilator \cite{Elowitz} is a synthetic genetic regulatory
network in Escherichia Coli consisting of three genes ($lacI, tetR, cI$) connected in a feedback loop such
that each gene inhibits the next gene and is inhibited by the previous
one. The repressilator exhibits a stable oscillation with fixed time
period that is witnessed by the synthesis of a green fluorescent protein ($GFP$).
%%
According to \cite{Elowitz}: \emph{``the resulting oscillation with typical periods of hours are slower than the cell division cycle so the state of the oscillator has to be transmitted from generation to generation''}.  

It may be modelled with five entities $\entities =\{lacI, tetR, cI, GFP, gen\}$.
The first three represent the three genes,  each with two levels denoting an \emph{on}/\emph{off} state where level
$1=on$ has decay time $2$. $GFP$ has two levels (\emph{on}/\emph{off}) with unbounded decay. 
Finally, the evolution of successive generations is represented by entity $gen$ with, let's say, seven levels (with unbounded decay) to differentiate among generations.  
 
The behaviour is modelled  with \potential and \obligatory activities: \potential activities (on the left below) are used to represent the feedback loop between $lacI, tetR,$ and $cI$. \Obligatory activities  (on the right) handle the synthesis of $GFP$ and show the evolution of generations.
%
\[
\begin{array}{ll}
 \alpha_1: \activ{\emptyset}{(lacI,1)}{2}{(tetR,+1)} \quad  &  \beta_1: \activ{(lacI,1)}{\emptyset}{0}{(GFP, +1)}  \\
 \alpha_2: \activ{\emptyset}{(tetR,1)}{2}{(cI,+1)} &  \beta_2: \activ{(lacI,0)}{(lacI,1)}{0}{(GFP,-1)}\\
 \alpha_3: \activ{\emptyset}{(cI,1)}{2}{(lacI,+1)} &  \beta_3: \activ{(gen,0)}{(gen, 6)}{1}{(gen, +1)}\\
 &  \beta_4: \activ{(gen,6)}{\emptyset}{1}{(gen, -6)} 
\end{array}
\]
In particular, the evolution of successive generations is circular modulo 7 (generations evolve from 0 to 6 and then they start again from 0) this allows us to have a finite number of levels for $gen$). Activities $\beta_3$ and $\beta_4$ describe such a behaviour. Notice in particular the use of inhibitor $(gen, 6)$ in $\beta_3$ that allows to apply the activity for all levels of $gen$ from 0 to 5 and to apply $\beta_4$ for level 6.


%
We now exhibit a possible execution scenario for the repressilator, the table below shows for each line the current level of each entity:
\[
\begin{tabular}{c|c|c|c|c|l}
 ~$lacI$~ & ~$tetR$~ & ~$cI$~ & ~$GFP$~ & $gen$ & ~description\\
   0    &  0     &    1 &   0   &   0   & ~initial state\\
   0    &  0     &    1 &   0   &   0   & ~after Ph. 2 (tick)\\
   0    &  0     &    1 &   0   &   0   & ~after Ph. 1 (no activity)\\
   0    &  0     &    1 &   0   &   1   &  ~after Ph. 2  (tick and activity $\beta_3$)  \\
   0    &  1     &   1  &   0   &   1   & ~after Ph. 1 (activity $\alpha_1$)\\
   0    &  1     &   1  &   0   &   1   & ~after Ph. 2 (tick)\\ 
   0    &  1     &   1  &   0   &   1   & ~after Ph. 1 (no activity)\\
   0    &  1     &   0  &   0   &   2   & ~after Ph. 2  (tick, decay of $cI$  and activity $\beta_3$)\\
   1    &  1     &   0  &   0   &   2   & ~after Ph. 1 (activity $\alpha_3$)\\    
   1    &  1     &   0  &   0   &   2   & ~after Ph. 2 (tick)\\ 
   1    &  1     &   0  &   0   &   2   & ~after Ph. 1 (no activity)\\ 
   1    &  0     &   0  &   1   &   3   & ~after Ph. 2  (tick, decay of $tetR$ and activities $\beta_1,\beta_3$)
\end{tabular}
\]

\hfill $\diamond$
\end{example}



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\subsection{\bf Petri Net formalisation.}\label{sec:petri} 

The semantics of \andy networks is formalised in terms of high-level Petri nets. 
We recall here the  notations together with some elements of their
semantics  \cite{DBLP:series/eatcs/Jensen92}.

\begin{definition}\label{def:PTS} 
A {\bf high-level Petri net}  is a tuple $(P, T, F, L, M_0)$ where:
\begin{itemize}
 \item $P$ is the set of places, $T$ is the set of transitions, with $P \cap T = \emptyset$;
 \item $F \subseteq (P \times T) \cup (T \times P)$ is the set of arcs; 
 \item $L$ is the labelling  of $P \cup T\cup F$ defined as follows:
\renewcommand{\labelitemii}{$-$} 
\begin{itemize}
 \item $\forall p \in P$, $L(p) $ is a set of values;% the type of $p$;%, \ie a (possibly infinite) set or Cartesian product of sets of integer values;
 \item $\forall t \in T$, $L(t)$ is a  computable Boolean expression (a guard);%  with variables and integer values; 
 \item and $\forall f\in F$, $L(f)$  is a  tuple of variables or values.% compatible with the adjacent place.  
\end{itemize}
 \item $M_0$ is the initial marking putting a multiset of tokens in $L(p)$ in each   $p \in P$.\looseness=-1
\end{itemize} 
\end{definition}
The behaviour of high-level Petri nets is defined as usual: markings (states) are functions from places in $P$ to multisets of possibly structured tokens in $L(p)$. A transition $t \in T$ is \emph{enabled} at  marking $M$, if there exists a valuation $\sigma$ of all variables in the labelling of $t$  such that the guard $L(t)$ evaluates to true ($L_{\sigma}(t) = \mathsf{true}$) and there are enough tokens in all input places $p$ of $t$ to satisfy the corresponding input arcs, \ie $L_{\sigma}((p,t)) \in M(p) $ for each input place $p$ of $t$.
% and all input arc labels are $ \qquad \forall q \in \pre{t},   [L((q,t))]_{\sigma} \in M(q),$$
%where $\pre{t} =\{q\in Q \mid (q,t) \in F\}$ denotes the set of input places of a transition $t$. 
The firing of $t$ produces the marking $M'$:  
$\forall p \in P, M'(p) = M(p) - L_{\sigma}((p,t)) + L_{\sigma}((t,p))$
with $L_{\sigma}(f) = 0 $ if $f\notin F$, $-$ and $ +$ are  multiset  operators for removal and adding of one element, respectively. We denote it by $\firing{M}{t}{\sigma}{M'}$. 
Observe that for \andy we are considering a subclass of high-level Petri nets where at most one arc per direction  for each pair place/transition   is allowed and only one token can flow through. So, starting from an initial marking associating one token per place, this situation is maintained at each step of the net evolution.

As usual, in figures, places are represented as rounds, transitions as squares, 
and arcs as directed arrows.
By convention, primed version of variables (\eg $x'$) are used to annotate 
output arcs of transitions, their evaluation is  
possibly computed using unprimed 
variables (\eg $x$) appearing on input arcs. 
With an abuse of notation, singleton markings are denoted without brackets, 
the same is used in  arc annotations. 
Also we refer to valuated variables without effectively mentioning the valuation $\sigma$: \eg we say that the current value of variable $w$ is $w$ instead of $\sigma(w)$.
An example of firing is shown in  Figure \ref{fig:firing}.


We say that a marking $M$ is \emph{reachable} from the initial marking 
$M_0$ if there exists a firing sequence $(t_1, \sigma_1), \mydots, (t_n, \sigma_n)$ such that 
 $\firing{M_0}{t_1}{\sigma_1}{M_1} \mydots \firing{M_{n-1}}{t_n}{\sigma_n}{M}$.
The state space of an initially marked Petri net 
comprises all the reachable markings.

\begin{figure}
\begin{center}
\begin{tikzpicture}[node distance=1.5cm,>=stealth',bend angle=45,auto]

  \node [place, label = left: $p_1$] at (-1.5,0) (Q1){$7$};
  \node [place, label=left:$p_2$]at (0,-1.5)  (Q2){$5$};
  \node [place, label = right: $p_3$] at (1.5,0) (Q3){};
  
  \node [transition] (t) [  label=above : {$\begin{array}{c}
                                                {\scriptstyle x>y ~\wedge x'=x+y}
                                               \end{array}$}] at (0,0) {$t$}
   edge [pre]    node[below] {$x$}    (Q1)
   edge [pre]    node[right] {$y$}    (Q2)
   edge [post]   node[below]{$x'$}    (Q3);
\end{tikzpicture}
\hspace*{2cm}
\begin{tikzpicture}[node distance=1.5cm,>=stealth',bend angle=45,auto]

  \node [place, label = left: $p_1$] at (-1.5,0) (Q1){};
  \node [place, label=left:$p_2$]at (0,-1.5)  (Q2){};
  \node [place, label = right: $p_3$] at (1.5,0) (Q3){12};
  
  \node [transition] (t) [  label=above:{ $\begin{array}{c}
                                                {\scriptstyle x>y ~\wedge x'=x+y}
                                               \end{array}$}] at (0,0) {$t$}
   edge [pre]    node[below] {$x$}    (Q1)
   edge [pre]    node[right] {$y$}    (Q2)
   edge [post]   node[below]{$x'$}    (Q3);
\end{tikzpicture}
\end{center}
\caption{\label{fig:firing}Example of firing of transition $t$ with $\sigma = \{ x= 7, y=5, x'=12\}$.}
\end{figure}

\paragraph{\bf \andy's formalisation.} In order to properly describe decay and the semantics of activities we need to record three kinds of information:
\begin{enumerate}
 \item (decay) since how long an entity is at the current level since the last update,
 \item (activators) since how long an entity is available at a level less or equal to the current one (recall that levels are inclusive),  
 \item (inhibitors) since how long a level greater than the current one has been left.
\end{enumerate}
Notice that duration relative to Items 1 and 2 for the current level can be different: \ie  if an entity is continuously sustained
by some activities it remains available in the system at the current
level for a period that may be longer than the corresponding decay
time. Thus  Item 1 reports the duration since the \emph{last} update (\ie when the activity concerning the entity has been performed) and Item 2 is the duration since the \emph{first} appearance of the entity at that level.

For this reason, we model each entity $\entity \in \entities$ by a single Petri net place $p_{\entity}$ carrying  tuples
$\tuple{\bool_\entity, \refr_\entity, \birth_\entity}$ as tokens where
%%
$\bool_\entity$ is  the current level of $\entity$; 
%%
$\refr_\entity$ is a counter storing the duration spent in the current level since the last update (Item 1): $0 \leq \refr_\entity \leq \life_\entity(\bool_\entity)$; 
%%
and $\birth_\entity$ is a tuple of counters with $\setlev_\entity$ fields (one for
each level).  Each counter $\birth_{\entity}[i]$ contains a
duration interpreted as follows:
$$
\birth_{\entity}[i] = 
\begin{cases}
\begin{minipage}{6.5cm}
  since how long $\entity$ has reached level  $i$ (Item 2)
\end{minipage}
 &  \text{for }  0 \leq i \leq \lev_{\entity}\\
\\
\begin{minipage}{6.5cm}
   since how long $\entity$ has left level $i$ (Item 3)\\(or $0$ if $\entity$ has not yet reached $i$).
\end{minipage} 
& \text{for } \lev_{\entity} < i \leq \setlev_{\entity}-1              
\end{cases}
$$

Each place is initialised to $\tuple{\level{\entity}, 0,
  0^{\setlev_\entity} }$ where $\level{\entity}$ is the given initial
level of $\entity$ and $0^{\setlev_\entity}$ denotes vector $\birth_\entity$ of counters
uniformly initialised with zero.
 
It is worth observing that once a duration in $\birth$ has reached the maximum duration for all activities ($D = \max\{\Delta_\rho \mid \rho \in \Syn \cup \Reac \}$) it is useless to increment it, as for all duration greater than $D$, the guards of all activities involving the entity are satisfied.  


\begin{newnotation}
      We denote
      by $\birth\sub{x}{[l \mydots l + n ]}$ the systematic update of
      values in counters $\birth[l] \mydots \birth[l+n]$ to $x$. 
      $\inc{D}{\birth}$ is the increment by one of all values in $\birth$:
      \ie $\forall k \in [0 \mydots \setlev_\entity-1]$, $\birth[k] =
      \min(\birth[k]+1, D)$ with $D \in \nat$.
\end{newnotation}

In order to cope with time aspects present in \andy, we introduce a tick transition $t_c$ (Figure \ref{fig:clock}) that represents the time progression.  Its effect is the one described by Phase Ph. \ref{item:clock} above: it increments  counters of entities,  takes care of the decay and simultaneously performs all enabled \obligatory activities. More precisely: % updates due to time progression  are performed as follows:
      \begin{description}
      \item [Time: ] All counters in the entities are incremented by
            one: we pass from $\tuple{\lev, \refr, \birth}$ to
            $\tuple{\lev, \refr+1, \inc{D}{\birth}}$.  %We assume that
            Counter $\refr$ remains unchanged for levels of unbounded
            duration.
      \item [Decay: ] An entity may stay at level $\lev$ for
            $\life(\lev)$ time units (\ie $\refr \leq \life(\lev)$).
            Decay happens as soon as the interval $\life(\lev)$ is
            elapsed and is obtained by decreasing the level by one, till
            reaching level zero.  %More formally, 
            We pass from
            $\tuple{\lev, \refr, \birth}$ to $\tuple{\lev - 1, 0,
              \inc{D}{\birth}\sub{0}{\lev}}$.

      \end{description}

Updates on a place $p_r$ of a result entity $r$ at level $\lev_r$ containing token  $\tuple{\lev_r, \refr_r, \birth_r}$, due to an activity $\rho$ (\ie $(r,\pm n) \in R_{\rho}$ with $n\in \nat$) 
are performed\footnote{For the sake of simplicity, in the explanation we consider updates that do not exceed allowed boundaries, these cases are handled as expected in Definition \ref{def:rs}.} as follows : 
\begin{description}
 \item[Increase $+n$:]  
   $\tuple{\lev_r + n, 0, \birth_r\sub{0}{[\lev_r + 1\mydots \lev_r + n]}}$.
  $\lev_r$ becomes $\lev_r +n$, $\refr_r$ is reset to 0 restarting the counter  for level $\lev_r +n$ and 
$\birth_r$ is updated to record the change to the new level by resetting the counters $\birth_r[i]$ for all levels $i$ from $\lev_r + 1$ to $\lev_r + n$.
 
 \item[Maintenance $0$:]  $\tuple{\lev_r, 0, \birth_r}$: the current level is recharged resetting $\refr_r$ to 0.

  \item[Decrease $-n$:]  
	    $\tuple{\lev_r - n, 0, \birth_r\sub{0}{[\lev_r-n+1 \mydots \lev_r] }}$. This case is symmetric to the increase one, except for the treatment of $\birth_r$ that is updated by resetting  counters   $\birth_r[i]$
for all levels $i$  from  $\lev_r - n + 1$ to $\lev_r$.
\end{description}
%%
%%
% Once performed, the next execution of the same activity is possible after at least $\dur{\rho}$ time units.






Finally, every \potential activity $\alpha \in \Reac$ is modelled with  a transition $t_{\alpha}$ (Figure \ref{fig:reaction})  %with a special place $p_{\alpha}$ that  ensures that the same activity is not executed more than once every $\dur{\alpha}$ time units, 
complying to Phase Ph. \ref{item:reaction} above.







%\subsection{\bf Petri net modeling of \andy}


\begin{figure}[t]
\centering
\subfigure[Clock transition. \label{fig:clock}]{
\begin{tikzpicture}[xscale=1.2, yscale=1.2, node distance=1.5cm,>=stealth',bend angle=20,auto]

  \node [place, label = right: $p_{\rho}$, label = above: $\vdots$,label = below: $\vdots$] at (3.5,0) (Pc){};
  \node [place, label=left:$p_{\entity}$,label = above: $\vdots$,label = below: $\vdots$]at (0,0)  (R){};

  \node [transition] (tC) at (2,0) {$t_{c}$}
   edge [pre,bend left]    node[above] {$w_{\rho}$}    (Pc)
   edge [post,bend right]   node[below]{$w'_{\rho}$}    (Pc)
   edge [pre,bend right]    node[above] {${\scriptstyle \tuple{\bool_{\entity},\refr_{\entity},\birth_{\entity}}}$}    (R)
   edge [post,bend left]    node[auto] {${\scriptstyle \tuple{\bool_{\entity}',\refr_{\entity}',\birth_{\entity}'}}$}    (R);
  
\end{tikzpicture}
} \quad \quad 
\subfigure[Transition for $\alpha\in\Reac$. \label{fig:reaction}]{
\begin{tikzpicture}[node distance=1.5cm,>=stealth',bend angle=20,auto]

 \node [place, label = right: $p_r$,label = left: $\cdots$] at (0,-1.5) (P){};
 %\node [place, label = above: $\Pclock$] at (0,2) (Pc){};
 \node [place, label = below: $p_a$,label = above: $\vdots$] at (-2,0) (R){};
 \node [place, label = below: $p_i$,label = above: $\vdots$] at (2,0) (I){};
 \node [place, label = right: $p_{\alpha}$] at (0,1.5) (qw){};
 
 \node [transition] at (0,0) (tr) {$t_{\alpha}$}
  % edge [pre and post]    node[above] {$z$}    (Pc)
   edge [pre and post]    node[above] {${\scriptstyle \tuple{\bool_{\entity},\refr_{\entity},\birth_{\entity}}}$}   (R)
   edge [pre and post]    node[above] {${\scriptstyle \tuple{\bool_{\entity},\refr_{\entity},\birth_{\entity}}}$}     (I)
   edge [pre, bend right]    node[left] {${\scriptstyle \tuple{\bool_r,\refr_r,\birth_r}}$}   (P)
   edge [post, bend left]   node[right] {${\scriptstyle \tuple{\lev'_r,\refr'_r,\birth_r'}}$}     (P)
   edge [pre, bend left]    node[left] {$w_{\alpha}$}   (qw)
   edge [post, bend right]   node[right] {$w'_{\alpha}$}     (qw);
   
  
\end{tikzpicture}
}
\caption{Scheme of Petri net modelling of \andy.}
\end{figure}



\begin{definition}\label{def:rs}
 Given an  \andy network $(\entities, \Syn, \Reac)$ with initial 
state $(\entity, \level{\entity})$ for each $\entity \in \entities$, 
%Let $D = \max\{\dur{\alpha} \mid \alpha \in \Syn \cup \Reac\}$. 
the high-level Petri net representation is defined as 
tuple $(P, T, F, L, M_0)$ where  $\lev$, $\lev'$, $\refr$, $\refr'$, $\birth$, 
$\birth'$, $w$, $w'$ range over variables and:
\begin{itemize} 
 \item $P = \{p_\entity \mid  \entity\in{\entities}\} \cup \{p_{\rho} \mid  \rho \in \Syn \cup \Reac \}; $
 \item $T= \{t_c \} \cup \{t_{\alpha} \mid \alpha \in \Reac \}$;

 \item $F= \{(p,t_c),(t_c, p) \mid p \in P \} \quad  \cup$ \\
 \hspace*{0.63cm} $ \{ (p_\entity,t_{\alpha}),(t_{\alpha}, p_\entity), (p_{\alpha},t_{\alpha}),
		(t_{\alpha}, p_{\alpha}) \mid \alpha \in \Reac, 
		\entity \in A_{\alpha} \cup I_{\alpha} \cup R_{\alpha} \}$

 \item Labels for places in $P$: 
\[
\begin{array}{lr}
  L(p_{\rho})  = [0 \mydots D] & \text{ for each } \rho \in \Syn \cup \Reac  \\
L(p_\entity) = [0\mydots\setlev_\entity-1] \times [0 \mydots d] \times  [0 \mydots D]^{\setlev_\entity}  
& \text{ for each } \entity\in\entities 
\end{array}
\]
with $d = \max\{\life_{\entity}(i)\mid i \in [0 \mydots \setlev_\entity-1]\}$ and $D = \max\{\dur{\rho}\mid \rho \in \Syn \cup \Reac\}$
\item Labels for arcs in $F$:
\[
\begin{array}{llr}
 L((p_{\rho}, t_c)) = w_{\rho} & L((t_c, p_{\rho} )) = w'_{\rho} & \text{ for each } \rho \in \Syn \cup \Reac\\
    L((p_\entity, t_c)) = \tuple{\bool_\entity, \refr_\entity, \birth_\entity} 
		& L((t_c, p_\entity)) = \tuple{\bool_\entity', \refr_\entity', \birth_\entity'} 
		& \text{ for each } \entity\in\entities
\end{array}
\]
For each \potential activity $\alpha\in\Reac$ and $\entity \in A_{\alpha} \cup I_{\alpha} \cup R_{\alpha}$:
\[
\begin{array}{llr}
      L((p_\entity, t_{\alpha})) = \tuple{\bool_\entity, \refr_\entity, \birth_\entity} \quad 
		& L((t_{\alpha}, p_\entity)) = \begin{cases}
                           \tuple{\bool_\entity, \refr_\entity, \birth_\entity} & \text{if } \entity\notin R_{\alpha}\\
                           \tuple{\bool_\entity', \refr_\entity', \birth_\entity'} & \text{otherwise}
                          \end{cases}\\

     L((p_{\alpha}, t_{\alpha})) = w_{\alpha} & L((t_{\alpha}, p_{\alpha} )) = w'_{\alpha}
\end{array}
\]

\input{algorithm}


\item Labels for transitions in $T$:
 \[
 L(t_c) = \mathit{Calc}
 \]
where $\mathit{Calc}$ is a complex formula given by Algorithm  \ref{fig:calc}  that contains  equalities of the form $x_j' = c_j$ for each net variable $x_j'$ labelling an output arc  of $t_c$.
The algorithm takes as input all possible valuations $\sigma$ for net variables of input arcs of $t_c$.
 
\[
\begin{array}{lcl}
L(t_{\alpha}) \! &=&  w_{\alpha} \geq \dur{\alpha}  \wedge  w'_{\alpha} = 0 \quad \wedge \\
&& \bigwedge_{(a, \level{a}) \in A_{\alpha}}(\lev_a \geq \level{a} \wedge \birth[\level{a}] \geq \dur{\alpha} )  \quad \wedge\\
&&  \bigwedge_{(i, \level{i}) \in I_{\alpha}}(\lev_i < \level{i} \wedge \birth[\level{i}] \geq \dur{\alpha} )  \quad  \wedge \\


&& \bigwedge_{(r,+n) \in R_{\alpha}} C_{r+} \wedge \bigwedge_{(r,0) \in E_{\alpha}}C_{r0} \wedge  \bigwedge_{(r,-n) \in R_{\alpha}}C_{r-},\\  
\end{array}
\]
for each \potential activity $\alpha\in\Reac$ with
\[
\begin{array}{ll}
C_{r+} \!&:\!  \tuple{\lev'_r, \refr'_r, \birth'_r} = \tuple{\min(\lev_r +n, \setlev_r-1), 0, \birth_r\sub{0}{[\lev_r +1 \mydots \min(\lev_r +n, \setlev_r-1 ]}} \\
C_{r0} \!&: \! \tuple{\lev'_r, \refr'_r, \birth'_r} = \tuple{\lev_r, 0, \birth_r} \\

C_{r-} \!&: \! \tuple{\lev'_r, \refr'_r, \birth'_r} = \tuple{\max(0,\lev_r-n), 0, \birth_r\sub{0}{[\min(0,\lev_r-n) \mydots \lev_r]}}.
\end{array}
\]

 

\item The initial marking is $M_0(p_{\entity})= \tuple{\level{s},0,0^{\setlev_\entity}}$ for each  $\entity\in \entities$ and $M_0(p_{\rho})= 0$ for 
$\rho\in \Syn \cup \Reac$.
% , the initial marking $M_0$ is:
% \[
% M_0(p) =
% \begin{cases}
%  0 & \text{if } p = p_{\rho} \\
%  \tuple{\level{s},0,0^{\setlev_\entity}}& \text{if } p = p_\entity
% \end{cases}
% \]
\end{itemize}
\end{definition}



We now comment on the transitions guards. 
Input and output arcs between the same place and transition with the same label (read arcs) are denoted in figures with a double-pointed arrow with a single label.
 We detail the conditions and the results of firing of $t_{\alpha}$ corresponding to a \potential activity $\alpha= \activ{A}{I}{\dur{}}{R}$. We have:
\begin{itemize}
\item place $p_{\alpha}$ retains whether the activity can be performed $w_{\alpha}\geq\dur{}$ and after firing the token is updated to  $0$, ($w'_{\alpha}=0$).
 \item each activator $a \in A$ has to be present at least at  level $\level{a}$ for at least $\dur{}$ time units, this is expressed by  guard $\lev_a \geq \level{a} \wedge \birth[\level{a}] \geq \dur{} $; 
 \item each inhibitor $i\in I$ has not to exceed  level $\level{i}$ for at least $\dur{}$ time units,  this is guaranteed by guard 
$\lev_i < \level{i} \wedge \birth[\level{i}] \geq \dur{}$;
\item each result is updated according to conditions $C_{r+}$ (increase), $C_{r0}$ (maintenance), $C_{r-}$ (decrease);

%
%if a result $r \in R $ has to be increased maintained or decreased then the corresponding place $p_r$ is updated to:
%\begin{description}
 %\item [increase:] $\tuple{\min(\lev_r +n, \setlev_r-1), 0, \birth_r\sub{0}{[\lev_r +1 .. \min(\lev_r +n, \setlev_r-1 ]}}$
 %\item [maintenance:] $ \tuple{\lev'_r, \refr'_r, \birth'_r} = \tuple{\lev_r, 0, \birth_r} $
 %\item [decrease:] $ \tuple{\lev'_r, \refr'_r, \birth'_r} = \tuple{\max(0,\lev_r-n), 0, \birth_r\sub{0}{[\min(0,\lev_r-n).. \lev_r]}}$. 
%\end{description}
\end{itemize}
Transition $t_c$ is always enabled, the effect of the firing is that  all connected places are modified according to the valuation $\sigma$ ``computed'' by algorithm $\mathit{Calc}$. The firing updates the levels of all entities dependently if they are involved in a \obligatory activity or decay and increments the counters in $p_{\rho}$ ($w'_{\rho} = w_{\rho}+1$) for all $\rho \in \Syn \cup \Reac$ unless a \obligatory activity $\beta$ is performed in which case the corresponding $w'_{\beta}$ is set to 0.\looseness =-1

\begin{figure}[t]
 \centering
 \includegraphics[width=\textwidth]{repress.png}
\caption{The general shape of the \andy network for the Repressilator from Example  \ref{ex:rep} (arcs annotations are omitted).}\label{fig:repress}
\end{figure}


\begin{example}\label{ex:rep-petrinet}

Figure \ref{fig:repress} gives the Petri net representation of the \andy system introduced in Example \ref{ex:rep},
obtained via the Snakes library.

\end{example}

\paragraph{\bf Expressiveness.}
The first result that we show for \andy networks concerns the finiteness of the state space. 
Actually, even if the number of species and levels is finite, as a state of the network includes an information concerning time (that could be unbounded, if one takes
a naive representation of dates), this could lead to an infinite state space. 
Nonetheless in our case, when the local counters have reached the maximal duration for activities $D$, 
no new behaviours can be triggered by the system. 
This immediately provides a way of abstracting without loosing precision and thus getting a symbolic representation. 

\begin{theorem}
 Let $N$ be an \andy network. The state space of $N$ is finite
\end{theorem}
\begin{proof}[Sketch]
 The proof follows  from the observation that the type of tokens for each place is finite and each place can have at most one token. \qed
\end{proof}

Following from this result, it is easy to build a timed automata model \cite{DBLP:journals/tcs/AlurD94} of the semantics of \andy, this could help in better understanding and supporting the notion of time. Since their introduction timed automata have been widely used to model complex real time systems and to check their temporal properties. Since then many variants have been considered depending of
synchronisation policies (synchronous, urgent, broadcast \cite{uppaal}).

The timed automaton should have as many localities as the markings of the Petri net. Tick transitions of the Petri net are represented by the natural time progression in timed automata and the discretization is imposed by specific guards. This simple encoding has the same size as the state space and therefore exponential. In Appendix \ref{sec:timedauto} we provide a  different encoding that is compositional with respect to entities and \potential activities. This encoding is more compact but still exponential because of the nature of \obligatory activities that have to be treated simultaneously and dynamically.  

Next theorem states that the semantics in terms of timed
automata is equivalent to the one of Definition \ref{def:rs}.

\begin{theorem}\label{th:correct}
\andy can be encoded into timed automata preserving the semantics. The size of the encoding is exponential.
\end{theorem}

This last result not only clarifies that \andy can suitably handle time aspects but also it underlines the concision of the model. 
As a matter of fact, the syntax of \andy profits from the advantages of Petri net modelling: \ie the state space is the combination of local information.  This is  not the case for automata where the number of localities is already exponential.
Moreover the rule-based architecture of \andy guarantees a modular construction of systems which is not the case neither (in general) for Petri nets nor for timed automata. 
     
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Application to toxicology}\label{sec:toxic}
In this section, we explore the notion of toxicity and provide a methodology to analyse it in \andy.
Toxicology~\cite{Waters2004} studies the adverse effects of the exposures to chemicals at various levels of living entities: organism, tissue, cell and intracellular molecular systems. During the last decade, the accumulation of genomic and post-genomic data together with the introduction of new technologies for gene analysis has opened the way to  \emph{toxicogenomics}. Toxicogenomics combines toxicology with
``Omics'' technologies\footnote{``Omics'' technologies are methodologies such as genomics, transcriptomics, proteomics and metabolomics.}
to study the \emph{mode-of-action} of toxicants or environmental stressors on biological systems. The mode-of-action is understood as the sequence of events from the absorption of chemicals to a toxic outcome.
Toxicogenomics potentially improves clinical diagnosis capabilities and facilitates the identification of potential toxicity  in drug discovery~\cite{Foster2007} or in the design of bio-synthetic entities~\cite{synthetic}. 

The main approach used in toxicogenomics employs empirical analysis like in the identification of molecular biomarkers, \ie indicators of disease or toxicity in the form of specific gene expression patterns~\cite{DeCristofaro2008}.
Clearly, biomarkers remain observational indicators linking genes related measures to toxic states. In this proposal, we complement these empirical methods with a computational method  that aims at discovering the molecular mechanisms of toxicity.
This  way, instead of studying the  phenomenology of the toxic impacts, we focus on the processes triggering adverse effects on organisms.
Usually, the toxicity process  is defined as a sequence of physiological events that causes the abnormal behaviour of a living organism  with respect to its healthy state.

Healthy physiological states generally correspond to homoeostasis, namely a process that maintains a dynamic stability of internal conditions against changes in the external environment.
Hence, we will consider toxicity outcomes as deregulation of homoeostasis processes, namely deviation of some intrinsic referential equilibrium  of the system. 
Biological processes are usually given in terms of pathways which are  causal chains of the responses to stimuli, this way the deregulation of homoeostasis  appears as the unexpected activation or inhibition of existing pathways. Moreover,
in the context of toxicogenomics it is crucial to take into account at least two other parameters: the exposure time  and the thresholds dosage delimiting the ranges of safe and hazardous effects.  
Considering the qualitative nature of \andy models, we believe that the potential contribution lies in the analysis of the etiology of toxic processes. From this perspective, toxicity properties can be classified into two categories:
\begin{enumerate}
     \item \label{symptom} properties characterising pathological states or more generally associated to symptoms  
and
    \item \label{pathway} properties characterising sequences of states (traces), \eg non viable behaviours.
\end{enumerate}
The former class of properties basically leads to check the reachability of some states, while the latter may be used to unravel sequences of events leading to toxic outcomes. 

\andy allows us to describe such scenarios and the Petri net representation of a reaction network with the associated marking graph may be used to detect and predict toxic behaviours related to the dynamics of biological networks, both from the point of view of the dose (\eg some species exceed a pathological level) and from the point of view of the exposure (\eg some species persist in time beyond some maximal exposure time). Remark that defining toxicity by characterising properties on sequences of events is more general than just characterising pathological states. It addresses the needs to characterise multifactorial scenario (\eg periodical behaviours relating possibly several species) and to discover complex etiology among reactants. 

\paragraph{\bf From healthy states to healthy behaviours.}
In \andy a state is a marking in the marking graph, it collects for each species its current level and the corresponding $\lambda$ vector. The first step is to observe that some states can be naturally characterised as dangerous, they highlight particular symptoms (e.g. the dosage of a poison has become relevant). Let \danger be the set of such states. Symmetrically to dangerous states, there are some states that characterise a ``safe'' condition (e.g. the organism is functioning normally and does not presents hazardous symptoms), denoted \healthy. Notice that there are states that are not assigned to any of the two classes, we group them into the \remain set.
Characterising harmful states corresponds only to verify toxicity problems as in Item \ref{symptom}  above, more complex properties (Item \ref{pathway}) need to refer to (sets of) traces in the marking graph. We therefore classify possible behaviours starting from an initial state into four toxicity scenarios: 

\begin{enumerate}
\item those that lead to  states in \danger, 
\item those that irremediably quit \healthy and never reach it again: the tail of the trace visits only states in \danger or \remain (\eg when a gene mutation knock out some  pathway);

\item those that quit \healthy for too long in terms of ticks (even if they may reach it again), the regulation system allows, in principle, to reestablish the homeostatic equilibrium but the process takes too long, compromising the viability of the organism (\eg DNA damage caused by hydrogen peroxide in Escherichia Coli can be recovered to some extent by DNA repair enzymes but the perturbation may lead to cell apoptosis if to important or too long \cite{Imlay03061988});

% \item those that quit \healthy and reach a state too dissimilar from the ones in \healthy (\eg 
% a mutation) this may require the introduction of a notion of distance in the marking graph 
% that may depend on the specific domain under analysis,

\item those that quit \healthy an unbounded number of times (systemic metabolic dysfunction, slow poisoning) 

\end{enumerate}

All these properties can be given in terms of temporal logic formulae. As mentioned above,  the first item corresponds to reachability analysis thus instead of logic formulae ad hoc algorithms can be used as well. The second item can be dealt with formulae like the following in CTL: EF EG $ (s \notin \healthy)$, which means it is possible (EF) that starting from a certain moment all visited states (EG) will not belong to set \healthy. The remaining items correspond to more complex formulae whose formulation depends on the specific property of the system under consideration, for instance taking into account a precise number of ticks or a concrete sequence of states. 

\paragraph{\bf From healthy behaviours to healthy states.}
Notice that sometimes it is not sufficient or possible to describe the healthy condition of an organism (i.e. the ability of maintaining its capabilities) only in terms of the sets \healthy and \danger. Indeed, for some organisms their viability conditions are characterised by complex behaviours such as specific traces or oscillations. We, thus, define toxic (and symmetrically viable) behaviours $\mathcal{T}_{toxic}$ (resp. $\mathcal{T}_{viable}$) by giving the property characterising them (or simply by enumerating the traces).  The properties can be described by specific temporal logic formulae (that can be rather complex employing also recursive operators). For instance, we can consider viable all traces that infinitely repeat an ordered sequence of three states  as in our example of the repressilator (\eg only lacL ON, then only tetR ON, and finally only cI ON)  and where the period of oscillation is six ticks. 
This way toxicity can be checked by testing whether a trace belongs or not to $\mathcal{T}_{toxic}$.  
For instance, suppose that we want to describe the periodic activation of three genes $a, b$ and $c$, we could use this $\mu$ calculus formula to describe $\mathcal{T}_{viable}$:
$$\nu.X (a \wedge \langle-\rangle\langle-\rangle\langle-\rangle tt \wedge [-](b \wedge [-](c \wedge [-]X )))$$
which means we are in a state with $a$ activated, followed by three transitions, for each first transition we should have $b$ activated and then $c$ and  so on recursively. 
In our case, as the state space is finite, any \andy network may be simulated by a B\"uchi automaton and the questions above may correspond either to check word acceptance or language inclusion, which are both decidable questions in our case.


Moreover once we have given the specification of viable and toxic behaviour, we may infer several interesting sets of states:
\begin{itemize}
\item  states starting from which all possible futures are behaviours in $\mathcal{T}_{viable/toxic}$ 
\item states that have at least one possible future in $\mathcal{T}_{viable/toxic}$,
\item states that are included only in $\mathcal{T}_{viable/toxic}$.
\end{itemize}








\begin{example}
We generalise here the example of the repressilator given above: we take three entities 
$\{A, B, C\}$ each with decay duration 2 that are inhibited in turn by potential activities:
\[
\begin{array}{l}
 \alpha_1: \activ{\emptyset}{(A,1)}{2}{(B,+1)} \\
 \alpha_2: \activ{\emptyset}{(B,1)}{2}{(C,+1)} \\
 \alpha_3: \activ{\emptyset}{(C,1)}{2}{(A,+1)}  
\end{array}
\]
we consider a new perturbed variant of the system with an additional entity $D$ with only one level (\ie the entity is permanent) such that $\alpha_1$ and $\alpha_3$ remains unchanged and $\alpha_2$ is modified in  $\alpha_2': \activ{(D,1)}{(B,1)}{4}{(C,+1)}$. 

 \begin{figure}[h]
 \centering
 \begin{tikzpicture}[>=latex',xscale=.4, yscale=.5]
 
 \draw[draw=none,fill=blue!10] (-1,1) rectangle (22,2) ; 
 \draw[draw=none,fill=blue!10] (-1,3) rectangle (22,4) ; 
 
 \node at (-1,1) (p0) {};
 \node[left] at (p0.west) {$0$};
 \node at (-1,2) (p1) {};
 \node[left] at (p1.west) {$1$};
 \node at (-1,3) (c0) {};
 \node[left] at (c0.west) {$0$};
 \node at (-1,4) (c1) {};
 \node[left] at (c1.west) {$1$};
 
  \node[label=left:{$A'=\left[ \phantom{\begin{array}{l} a\\ a \end{array}} \right.$}] at (-1,3.5) (c) {};
 \node[label=left:{$A=\left[ \phantom{\begin{array}{l} a\\ a \end{array}} \right.$}] at (-1,1.5) (p) {};
 
 %\draw[->] (p0) -- (p1) ;
 %\draw[->] (c0) -- (c1) ;
 %\draw[->] (g0) -- (g1) ;
 %\draw[->] (b0) -- (b3) ;
 
 %points horizontaux
 \node at (-2,0) (t0) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:0] at (0,0) (t1) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:2] at (2,0) (t2) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:4] at (4,0) (t3) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:6] at (6,0) (t4) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:8] at (8,0) (t5) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:10] at (10,0) (t6) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:11] at (12,0) (t7) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:12] at (14,0) (t8) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:13] at (16,0) (t9) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:14] at (18,0) (t10) {};
 \node[rectangle,draw,fill=blue,scale=.3,label=below:16] at (20,0) (t11) {};
 
 \node[label=below right:$time$] at (22,0) (t12) {};
 
 % points verticaux
 \node[rectangle,draw,fill=blue,scale=.3] at (-1,1) (k1) {};
 \node[rectangle,draw,fill=blue,scale=.3] at (-1,2) (k2) {};
 \node[rectangle,draw,fill=blue,scale=.3] at (-1,3) (k3) {};
 \node[rectangle,draw,fill=blue,scale=.3] at (-1,4) (k4) {};
 
 \draw[->] (t0) -- (t12) ; %axe du temps
 \draw[->] (-1,0) -- (-1,5.5) ; %axe vertical
 
 \path[draw,thick] (0,2) --(2,2) --node[midway,above,sloped]{} (2,1) -- (6,1) --node[midway,above,sloped]{} (6,2) -- (8,2)  --node[midway,above,sloped]{} (8,1) -- (12,1) -- (12,2) -- (14,2) -- (14,1) -- (18,1) -- (18,2) -- (20,2) -- (20,1) -- (21,1);
 \node[right] at (21,1) {$\cdots$};
 \path[draw,thick] (0,4)--(2,4)  --node[midway,above,sloped]{} (2,3) -- (8,3)  --node[midway,above,sloped]{} (8,4) -- (10,4) -- (10,3) -- (16,3) -- (16,4) -- (18,4) -- (18,3) -- (21,3);
 \node[right] at (21,3) {$\cdots$};
 
 \end{tikzpicture}
 \caption{Entities levels evolving with time corresponding to the given scenario.}
 \label{fig:diagram}
 \end{figure}


Figure \ref{fig:diagram} above shows the periodic activation of entity $A$ (first line) each six time units while the second line shows the periodic behaviour of $A$ (second line, here denoted with $A'$ to avoid confusion) in the system perturbed with entity $D$. 
 

We now define $\mathcal{T}_{viable}$ as all traces where entity $A$ is activated every six ticks, clearly with this property all   traces in the  perturbed system will be toxic.

\hfill$\diamond$
\end{example}

Notice that formulae describing properties can be rather long requiring, for instance, a complex combination of all possible cases and interleaving. To this aim, in the future we plan to design a language with dedicated primitives that abbreviate common use cases (\ie macro instructions).
A second example coming from a different domain and involving a greater number of activities and entities  is given in Appendix \ref{sec:example}. 







% 
% The former class of properties basically consists in verifying reachability of some states, 
% while the latter  concerns pathways that highlight sequences of events leading to toxic outcomes.
% For instance, in the case of glucose regulation, we could verify whether glycemia levels are kept stable and whether they change in  case of ingestion of aspartame. More precisely, we could examine the causes and the symptoms of the hypoglycemia induced by the assimilation of aspartame. Hence  hypoglycemia is treated as a toxic state. 
% 


\section{Related works}\label{sec:related}

From a technical point of view, the closest related work is on reaction
systems \cite{DBLP:journals/ijfcs/BrijderEMR11} or their Petri net
representation \cite{koutnyreaction}. Although we use a similar
definition for activity, the semantics we have proposed is inherently
different: in \cite{DBLP:journals/ijfcs/BrijderEMR11} all enabled
reactions occur in one step while we consider two forms of activity evolution: simultaneous  for \obligatory activities and
interleaved for \potential ones.
%%
Furthermore, we have introduced discrete abstract  levels that, to the
best of our knowledge, are not taken into account in reaction systems. 
Also, we have presented a more elaborated notion of time that governs
both entities and activities. 
%%
In  \cite{DBLP:conf/birthday/BrijderER11} the authors consider an
extension of reaction systems with duration but it concerns only decay
and not  duration of activities. Furthermore the decay is  referred to
the number of steps and not to the discrete time progression as in our case.\looseness=-1  

From the Petri net modelling point of view, our representation of time is considerably different from the approaches traditionally used in time and timed Petri nets 
(\cite{DBLP:journals/tcs/CeroneM99} presents a survey with insightful
comparison of the different approaches).
%%
The main difference lies on the fact that the  progression of time is implicit and  external to the system. By contrast, in our proposal we have assumed the presence of  an explicit way of incrementing duration (modelled by synchronised counters). 
This is also  different from the notion of timestamps introduced in  \cite{timestamp} that again refers to an implicit notion of time.
Indeed,  our approach is conceptually closer to Petri nets with causal time  \cite{DBLP:journals/entcs/ThanhKP02} for the presence of an explicit transition for time progression. Nevertheless, in our approach time cannot be suspended under the influence of the environment (as is the case in \cite{DBLP:journals/entcs/ThanhKP02}).\looseness=-1
%%
 
In a broader sense, our work could also be related to P-systems
\cite{DBLP:journals/ijfcs/PaunPRS11,DBLP:journals/fuin/KleijnK11} or the
$\kappa$-calculus \cite{DBLP:journals/tcs/DanosL04} that  describe the
evolution of cells through rules. 
%%
Both these approaches are mainly oriented to simulation while we are
interested in  verification aspects. Finally, always related to the
modelling in Petri nets but with a different aim, levels have been used
in qualitative approaches to address problems related to the
identification of steady states in genetic networks such as in
\cite{DBLP:journals/nc/ChaouiyaNRT11}. 
%Similarly in workflow modeling, several models have been proposed \cite{DBLP:conf/dfg/EshuisW03}.
Nevertheless these contributions
abstract away from time related aspects that are instead central in our
proposal.  
Finally, it is interesting to notice that our formalism adds a new node in between timed and qualitative regions in the classification in   \cite{DBLP:conf/apn/HeinerG11}.



 
 
\section{Final remarks}\label{sec:concl}


We have introduced \andy %, a high-level Petri net framework
to model
time-dependent activity-driven systems, which  consist of a set of entities
present in the environment at a given level. Entities can age as time
passes by and their level is governed by a set of \obligatory and \potential activities with duration. 
The resulting Petri net is compact and provides an intelligible
graphical representation of the described activities. Moreover, the
network is easily scalable: the size of the network grows linearly with
the number of added entities and activities. 
%In an activity, entities can have the role of activators, inhibitors and
%results. 
%An activity has a given duration. 
 %and it can take place only if
%all activators are available and all inhibitors are not for the whole
%duration of the activity. Depending on the type of activity, results
%levels are either increased, maintained or decreased. 
We have shown that, despite the ability of modelling timed (thus infinite) systems, \andy networks have finite state space and that they have an alternative but exponential formulation with equivalent semantics into timed automata (TA).

The development of an equivalent semantics in terms of modular TA, was instrumental. Even if the TA formulation will have little impact on the final users, it helped us to clarify the semantics of the temporal relationships in \andy by confronting the high level Petri net approach with the direct representation of clocks accessible in TA. We insisted on a compositional form of TA, despite the fact that it will make the translation more complex than a direct translation, because we strongly advocate for modularity to face the scalability and reusability challenges posed by systems biology.

Finally we have discussed how toxicity problems can be addressed using our formalism. To the 
best of our knowledge this is the first attempt to give a general methodology, a systemic approach of toxicity analysis rather than a particular technique to verify them. 


As the semantics of \andy is given in terms of high-level Petri nets, \andy networks can be easily used as overlay of existing implemented tools.
%
%
%The main application of our work concerns the verification of properties
%of systems defined in terms of  activities. 
We have prototyped a
tool to simulate and build the state space of an \andy network 
%relying on the high-level Petri net semantics and 
using the SNAKES
toolkit~\cite{Pom-PNN-2008}. We
have also developed an alternative implementation using the MGS
framework~\cite{giavitto2008} to test several strategies that can be
used to avoid the combinatorial explosion and to build lazily the state
space. 
The first results are very positive,  \andy is particularly suited to describe biological systems, in particular regulatory networks and pathways whose formalisation is based on rules (activities). 
The tools above have been effectively used in the analysis of several
applications: in physiology (modelling the metabolism of glucose)  \cite{CDGHKFD14}, in the
modelling of regulatory networks (\eg the repressilator),
etc. %The results will be reported elsewhere.
We have used SNAKES and MGS as they allow a rapid prototyping but we plan to re-implement the framework and to base further developments on mainstream tools that  are well supported in the Petri net system biology community as simulation tool Snoopy \cite{DBLP:conf/apn/HeinerHLRS12} and related analysis tool MARCIE \cite{DBLP:conf/apn/HeinerRS13}.


%
 We plan to develop the \andy
approach along several directions. We want to extend decay function to arbitrary maps in $\setlev$. The
problem is to generalise the notion of increase and decrease in a
relevant way for the application. More structure can also be given to
the notion of entity. Characterising an entity by several attributes
(instead of only one level) raises the question of the causal
relationships between these attributes and their update
strategy. Another aspect that will be addressed in the future is to
consider dynamic models where entities are created and deleted, as
 in~\cite{Giavitto2012}. 
%Finally, we are also working on more efficient
% translation towards timed automata to make possible the use of
% model-checking tools like UPPAAL \cite{uppaal} on more significant examples.
Finally, it would certainly be interesting to find a systematic way for approximating stochastic rates into activities duration and decays.



%% 
%Related: reaction systems, petri nets with causal time and petri nets with timestamps
%


 \paragraph{\bf Acknowledgements.}
 This work was supported by the French project ANR BLANC 0307 01 - SYNBIOTIC. We thank the anonymous reviewers of PNSE for their comments on the preliminary version of this paper. We are grateful to M. Heiner,  F. Pommereau and Y.-S. Le Cornec for useful suggestions and discussions. 
 %  

\bibliographystyle{abbrv}
\bibliography{biblio}

\newpage

\appendix
\input{appendix.tex}
\end{document}







\section{Repressilator}\label{app:repr}

Another example that is easily expressible with \andy is the
repressilator \cite{Elowitz}.  It is a synthetic genetic regulatory
network consisting of three genes connected in a feedback loop, such
that each gene inhibits the next gene and is inhibited by the previous
one. The repressilator exhibits a stable oscillation with fixed time
periods.
%%
We model it in \andy with three entities $\entities =\{s_1, s_2, s_3\}$
each with two levels representing an \emph{on}/\emph{off} state. Level
$1=on$ for all entities has decay time $\life_{s_i}(1)=1$; the set of
\potential activities is (no \obligatory activity are required):
%
\[
\begin{array}{c}
 r_1: \activ{\emptyset}{s_1}{1}{s_2}\\
 r_2: \activ{\emptyset}{s_2}{1}{s_3}\\
 r_3: \activ{\emptyset}{s_3}{1}{s_1}\\
\end{array}
\]
%
The SNAKES rendering of the repressilator is given in Figure~\ref{fig:repress}.

\begin{figure}
{\centering
\includegraphics[width=\textwidth]{repress.png}}
\caption{Petri net representation of the repressilator with Snakes}\label{fig:repress}
\end{figure}

Figure~\ref{fig:mgsB} shows the complete state graph where we identify
the characteristic oscillating behaviour of the repressilator
(represented by cycles). However, notice that there are also paths
that leave the cycles to return into state $3$ where the 3 entities are
not expressed (level $0$). 

\begin{figure}[ht]
\begin{center}
\includegraphics[width=.5\textwidth]{B.pdf}
\end{center}
\caption{State graph of the repressilator}\label{fig:mgsB}
\end{figure}

In order to have a better understanding of these unexpected behaviours,
we use the MGS tool to generate the alternative state graph given
Figure~\ref{fig:mgsA} representing the time progression. In fact, in
this graph, each arc corresponds to a clock tick (phase
Ph.~\ref{item:clock}) followed by a composition of \potential
activities (phase Ph.~\ref{item:reaction}).
%%
The periodic behaviour is clearly restricted to a subgraph (bold
cluster on Figure~\ref{fig:mgsA}); the other cluster (dashed on
Figure~\ref{fig:mgsA}) gathers the unexpected states.
%%
After analysing this graph, transitions going from the periodic
cluster to the other corresponds to evolutions with only time
progression and decay without any activity. It may be interpreted as
those cases that are unlikely to happen but still (at least in
principle) plausible which may explain the lack of robustness of such
a synthetic genetic regulatory network when  implemented in real scenarios.

\begin{figure}[ht]
\includegraphics[width=\textwidth]{A.pdf}
\caption{Time progression in the repressilator state space}\label{fig:mgsA}
\end{figure}

%% Moreover we have studied its state space with MGS obtaining a graph
%% that shows the periodic behavior typical of the repressilator. Figure
%% \ref{fig:mgsB} shows an extract of the complete state space where we
%% identify the oscillation (represented by cycles). In the complete
%% graph (Figure \ref{fig:mgsA}) notice that there are also paths that
%% leave the cycles. This corresponds to evolutions with only time
%% progression and decay without any activity and it may be interpreted
%% as those cases that are unlikely to happen but still (at least in
%% principle) plausible.  Thus the model is suitable to address
%% robustness matters as we can test the model against the appearance of
%% even the less likely to occur events.

%% \begin{figure}[ht]
%% \includegraphics[width=\textwidth]{B.pdf}
%% \caption{Extract of the state space of the repressilator}\label{fig:mgsB}
%% \end{figure}



\end{document}



\comment{%
        Initialement on disait : ``chaque entité possède un ensemble
        d'attributs'', ce qui n'est pas le cas vu la formalisation. Si
        une entité possède un seul attribut, alors il vaut peut-être
        mieux parler de \emph{variable} ou d'\emph{observable} que
        d'entité (car on parle des variables ou des observables d'un
        système et qu'habituelement une entité/objet est caractérisé par
        plusieurs attributs): ``a measurable property of a physical
        system''. %
      }



\comment{%%
  Why restraining to decay? For instance, with this formalism,
  I cannot fill in a bathtub. 

  Using an arbitrary $d_\entity$ for the evolution of the state of
  $\entity$ subsumes both increase and decrease. Note also that it
  simplifies greatly the modeling of a \emph{clock}: it is enough to
  consider a $d_\entity$ with a limit cycle. Note that this does not
  prevent to check that $d_\entity$ is monotonic (that is, induce a
  strict order on $\setlev$ when $x \leq d_\entity(x)$ for all $x \in
  \setlev$) and to use this property in the tools (\eg to cut
  the state space).  %%
} %%



\comment{%%
  Again: we can have states that are not levels. It will extend at no
  cost the applicability of the formalism. The decreasing function on
  levels can be modeled by a map $d_\entity : \setlev \mapsto \setlev$. We can
  insists to have an infinite delay for the fixed points $x$ of $d_\entity :
  \life_\entity(x) = \omega \Leftrightarrow d_\entity(x) = x$ and that
  $\level{\entity}$ is a fixed point. 
  %% 
}



 \comment{%%
  The purpose/use of $v$ must be given here. %%
  To generalize decay, replace $v \in \mathbb{Z}$ by $d : \setlev_\entity
  \mapsto \setlev_\entity$ and the addition/subtraction by the
  application of $d$. Think to $d$ as the successor/predecessor
  function. Requiring $d$ \emph{injective} to have $d^{-1}$ is to strong
  because of the fixed points. What is the weakest relevant
  assumption?  %%
}


\comment{%%
  Why is the discharge not related to the duration of the call? How to
  model it?
        
  If the phone is plugged, it does not discharge when there is a phone
  call.  %%
} 
